/**
 * 
 */
package edu.cmu.mism.dgjava.algorithm.impl;

import java.util.Random;

import org.eclipse.core.runtime.IProgressMonitor;

import edu.cmu.mism.dgjava.data.models.option.BaseType;
import edu.cmu.mism.dgjava.data.models.option.Equity;
import edu.cmu.mism.dgjava.data.models.option.MonteCarloSimulation;
import edu.cmu.mism.dgjava.data.models.option.OptionCalculationResult;
import edu.cmu.mism.dgjava.data.models.option.impl.OptionModelFactoryImpl;

/**
 * @author Ankit
 * 
 */
public class MonteCarloSimulationAlg extends OptionsAlg {

	private Equity underlyingModel;
	private MonteCarloSimulation pricingModel;
	private int sampleSize;
	private static Random randomNumberGenerator = new Random();

	private double[] normalDistributions;
	private double[] stockPrices;
	private double[] optionPrices;
	private double meanOptionPrice = 0;
	private double optionValue;
	
	private double actualOptionValue = 0.0;

	/**
	 * @param underlyingModel
	 * @param pricingModel
	 */
	public MonteCarloSimulationAlg(Object underlyingModel, Object pricingModel) {
		super(underlyingModel, pricingModel);
		// TODO Auto-generated constructor stub
	}

	/* (non-Javadoc)
	 * @see edu.cmu.mism.dgjava.algorithm.IOptionAlgorithm#calculate(org.eclipse.core.runtime.IProgressMonitor)
	 */
	@Override
	public OptionCalculationResult calculate(IProgressMonitor monitor) {
		OptionCalculationResult result = OptionModelFactoryImpl.eINSTANCE
				.createOptionCalculationResult();

		initializeOptionsVariable();
		calculateOption();
		actualOptionValue = optionValue;
		calculateGreekLetters();
		
		result.setPrice(actualOptionValue);
		result.setDelta(delta);
		result.setGamma(gamma);
		result.setRho(rho);
		result.setTheta(theta);
		result.setVega(vega);

		return result;

	}

	/**
	 * 
	 */
	private void initializeOptionsVariable() {

		underlyingModel = (Equity) getUnderlyingTypeModel();

		pricingModel = (MonteCarloSimulation) getPricingModel();

		stockPrice = underlyingModel.getStockPrice();
		riskFreeRate = underlyingModel.getRiskFreeRate() / 100.0;
		volatility = underlyingModel.getVolatility() / 100.0;

		strikePrice = pricingModel.getExercisePrice();
		timeToExpiration = pricingModel.getTimeToExercise();
		sampleSize = pricingModel.getSampleSize();

		type = pricingModel.getType();

		normalDistributions = new double[sampleSize];
		stockPrices = new double[sampleSize];
		optionPrices = new double[sampleSize];

	}

	/* (non-Javadoc)
	 * @see edu.cmu.mism.dgjava.algorithm.impl.OptionsAlg#calculateDelta()
	 */
	@Override
	protected void calculateDelta() {

		initializeOptionsVariable();
		stockPrice += 1;
		calculateOption();
		delta = (optionValue - actualOptionValue)/1.0;
	}

	/* (non-Javadoc)
	 * @see edu.cmu.mism.dgjava.algorithm.impl.OptionsAlg#calculateGamma()
	 */
	@Override
	protected void calculateGamma() {
		
		initializeOptionsVariable();
		stockPrice += 2;
		calculateOption();
		gamma = ((optionValue - actualOptionValue)/2 - delta)/1.0;

	}

	/* (non-Javadoc)
	 * @see edu.cmu.mism.dgjava.algorithm.impl.OptionsAlg#calculateRho()
	 */
	@Override
	protected void calculateRho() {

		initializeOptionsVariable();
		riskFreeRate += 0.01;
		calculateOption();
		rho = (optionValue - actualOptionValue);
	}

	/* (non-Javadoc)
	 * @see edu.cmu.mism.dgjava.algorithm.impl.OptionsAlg#calculateTheta()
	 */
	@Override
	protected void calculateTheta() {
		// TODO Auto-generated method stub
		theta = (stockPrices[2] - stockPrices[0])/2.0;
		theta /= 365.0;

	}

	/* (non-Javadoc)
	 * @see edu.cmu.mism.dgjava.algorithm.impl.OptionsAlg#calculateVega()
	 */
	@Override
	protected void calculateVega() {

		initializeOptionsVariable();
		volatility += 0.01;
		calculateOption();
		vega = (optionValue - actualOptionValue)/(0.01 * 100);
	}

	/* (non-Javadoc)
	 * @see edu.cmu.mism.dgjava.algorithm.impl.OptionsAlg#calculateOption()
	 */
	@Override
	protected void calculateOption() {

		for (int i = 0; i < sampleSize; i++) {

			normalDistributions[i] = calculateNormalDistribution();
			stockPrices[i] = calculateStockPrices(normalDistributions[i]);
			optionPrices[i] = calculateOptionPrices(stockPrices[i]);
			meanOptionPrice += optionPrices[i];
		}
		meanOptionPrice /= sampleSize;

		optionValue = meanOptionPrice
				* (Math.pow(Math.E, -(riskFreeRate * timeToExpiration)));

	}

	/**
	 * @param d
	 * @return
	 */
	private double calculateOptionPrices(double d) {
		// TODO Auto-generated method stub

		if (type.getValue() == BaseType.CALL_VALUE) {
			return (d - strikePrice) > 0 ? (d - strikePrice) : 0;
		} else {
			return (strikePrice - d) > 0 ? (strikePrice - d) : 0;
		}

	}

	/**
	 * @return
	 */
	private double calculateNormalDistribution() {

		return randomNumberGenerator.nextGaussian();
	}

	/**
	 * @param d
	 * @return
	 */
	private double calculateStockPrices(double d) {
		// TODO Auto-generated method stub
		double power = (riskFreeRate - (volatility * volatility) / 2)
				* timeToExpiration + volatility * d
				* Math.sqrt(timeToExpiration);
		double price = stockPrice * (Math.pow(Math.E, power));
		return price;
	}
	
	
	/*
	 * The last algorithm are used to calculated the cumulative normal distribution with mean 0 and SD as 1 and inverse of that.
	 * These algorithms can be used in Black Scholes Algorithm 
	 * 
	 * Source: http://www1.fpl.fs.fed.us/CDF_Normal.java
	 */

//	/**
//	 * @param number
//	 * @return
//	 */
//	private double cumulativeNormalDistribution(double number) {
//
//		double probability = 0.0;
//		double absoluteValue;
//		double exponential;
//		double pdf;
//
//		final double p0 = 220.2068679123761;
//		final double p1 = 221.2135961699311;
//		final double p2 = 112.0792914978709;
//		final double p3 = 33.91286607838300;
//		final double p4 = 6.373962203531650;
//		final double p5 = .7003830644436881;
//		final double p6 = .3526249659989109E-01;
//
//		final double q0 = 440.4137358247522;
//		final double q1 = 793.8265125199484;
//		final double q2 = 637.3336333788311;
//		final double q3 = 296.5642487796737;
//		final double q4 = 86.78073220294608;
//		final double q5 = 16.06417757920695;
//		final double q6 = 1.755667163182642;
//		final double q7 = .8838834764831844E-1;
//
//		final double cutoff = 7.071;
//		final double root2pi = 2.506628274631001;
//		
//		absoluteValue = Math.abs(number);
//		if (number > 37.0) {
//
//	         probability = 1.0;
//
//	         return probability;
//
//	      }
//
//	      if (number < -37.0) {
//
//	         probability = 0.0;
//
//	         return probability;
//
//	      }
//	      
//	      exponential = Math.exp(-.5*absoluteValue*absoluteValue);
//
//	      pdf = exponential/root2pi;
//	      
//	      if (absoluteValue < cutoff) {
//
//	          probability = exponential*((((((p6*absoluteValue + p5)*absoluteValue + p4)*absoluteValue + p3)*absoluteValue +
//	              p2)*absoluteValue + p1)*absoluteValue + p0)/(((((((q7*absoluteValue + q6)*absoluteValue +
//	              q5)*absoluteValue + q4)*absoluteValue + q3)*absoluteValue + q2)*absoluteValue + q1)*absoluteValue +
//	              q0);
//
//	       } else {
//
//	          probability = pdf/(absoluteValue + 1.0/(absoluteValue + 2.0/(absoluteValue + 3.0/(absoluteValue + 4.0/
//	              (absoluteValue + 0.65)))));
//
//	       }
//
//	       if (number < 0.0) {
//
//	          return probability;
//
//	       } else {
//
//	          probability = 1.0 - probability;
//
//	          return probability;
//
//	       }
//
//	}
//	
//	/**
//	 * @param p
//	 * @return
//	 */
//	private double cumulativeNormalDistributionInverse(double p) {
//		
//		double arg,t,t2,t3,xnum,xden,qinvp,x,pc;
//
//	      final double c[] = {2.515517, 
//	      .802853,
//	      .010328};
//
//	      final double d[] = {1.432788,
//	      .189269,
//	      .001308};
//
//	      if (p <= .5) {
//
//	         arg = -2.0*Math.log(p);
//	         t = Math.sqrt(arg);
//	         t2 = t*t;
//	         t3 = t2*t;
//
//	         xnum = c[0] + c[1]*t + c[2]*t2;
//	         xden = 1.0 + d[0]*t + d[1]*t2 + d[2]*t3;
//	         qinvp = t - xnum/xden;
//	         x = -qinvp;
//
//	         return x;
//
//	      }
//
//	      else {
//
//	         pc = 1.0 - p;
//	         arg = -2.0*Math.log(pc);
//	         t = Math.sqrt(arg);
//	         t2 = t*t;
//	         t3 = t2*t;
//
//	         xnum = c[0] + c[1]*t + c[2]*t2;
//	         xden = 1.0 + d[0]*t + d[1]*t2 + d[2]*t3;
//	         x = t - xnum/xden;
//
//	         return x;
//
//	      }
//
//	}
	
	

}
